This invention relates to electrical capacitance tomography. Specifically, the invention relates to a new image reconstruction technique for imaging two- and three-phase flows using electrical capacitance tomography (ECT).
Applications of process tomography as a robust non-invasive tool for direct analysis of the internal characteristics of process plants in order to improve the design and operation of the equipment have increased in number in recent years (Williams and Beck, 1995, Beck and Williams, 1996). Process tomography involves utilization of tomographic imaging methods to manipulate data from remote sensors to obtain precise quantitative information from inaccessible locations. The need for tomography in industrial applications is analogous to the medical need for body scanners. A tomographic technique involves two basic tasks: (1) the acquisition of measurement signals from sensors located on the periphery of an object, such as a process vessel or column; and (2) the process of abstracting the measurement signals, which reveals information on the nature and distribution of components within the sensing zone, to form a cross-sectional image of the object. The task of generating image from the measurement signals is also known as tomographic reconstruction. The basic components of a tomographic instrument are, therefore, embodied in a sensor system, a signal/data acquisition system and a computer system for measurement control, image reconstruction and display.
Successful implementation of a tomographic technique lies on the selection of a sensor system deployed for the specified application and the tomographic image reconstruction algorithm suited to the sensor selected. A variety of sensing methods can be employed based on measurements of transmission, diffraction or electrical phenomena using radiation, acoustic or electrical sensors (Beck, 1995). For the purpose of multiphase flow imaging in process industries, one of the most critical things may be the speed of the technique to capture the real time data of the turbulent fluctuation in the flow field. In this regard, electrical capacitance tomography (ECT) is considered to be the most powerful tool among other available tomographic techniques because of its high-speed capability, low construction cost, high safety and suitability for small or large vessels. Currently commercially available electrical capacitance tomography systems could capture up to 100 image frames per second, compared to corresponding electrical resistance tomography (ERT) of 2 image frames per second, or X-ray computer tomography (CT) of only one image every 2 seconds.
ECT is gaining acceptance as a laboratory and industrial tool to analyze multiphase systems. Early efforts with ECT have been concerned with the imaging of two-phase stratified flows in industrial pipelines, especially oil-gas and oil-water flows from offshore oil wells (Yang et al., 1996), manufacturing processes involving gas-solid systems such as pneumatic conveyers (Ostrowski et al., 1997, Dyakowski et al., 1999) and gas-solid fluidizations (Dyakowski et al., 1997, Halow and Nicoletti, 1992, Wang et al., 1995), and trickle bed reactors for measuring water content (Reinecke and Mewes, 1997, 1998). However, application of the technique to more complex multiphase flow systems such as bubbly flows in gas-liquid as well as gas-liquid-solid systems which are widely used in chemical processes (Fan, 1989) is very limited. ECT has prospective uses for applications to gas-liquid and gas-liquid-solid systems in real chemical processes, since they are mostly using organic liquids (Fan et al., 1999), which are non-conductive, rather than water which is a widely used model liquid in laboratories. The implementation of the technique would considerably further research in these fields since the information provided by the rapid online imaging method offers a means to address long-standing problems in the modeling, optimization and control of these processes.
However, the most critical problems that still challenge the application of the technique to such systems may be the relatively low spatial resolution and the accuracy of the reconstructed image using existing techniques. Most work currently done by researchers is focused on efforts to address these problems. While leaving the problem on the spatial resolution to other researchers working on sensor hardware, this work deals with a reconstruction technique aimed at improving the accuracy of the reconstructed image. The selection of image reconstruction technique suited to the sensor is highly important, because it determines the quality of the image that gives the information required to analyze the flow system. In the case of capacitance tomography, unfortunately, the reconstruction problem is non-linear, so that commonly used and commercially available and well-developed reconstruction techniques for linear tomographies based on electromagnetic radiation as widely used in the medical field are not directly applicable to the non-linear problem. Although many studies have been reported on the development of image reconstruction techniques for capacitance tomography, the reconstructed results using the techniques reported so far are still more qualitative rather than quantitative. For application to relatively complex multiphase system such as bubbly flows, high accuracy of the reconstruction results is especially necessary. There is an urgent need for development of an accurate reconstruction technique for capacitance tomography to meet the requirement of real chemical process applications.
Multi-modal Tomographic Techniques for Three-phase Flow Imaging
Tomographic technique for two-phase flow imaging is referred to as single modal tomography, where the requisite sensed signal contains only one parameter in the object space, (e.g., energy absorption (electromagnetic radiations), permittivity (capacitance sensor) or conductivity (impedance sensor) (Nooralahiyan and Hoyle, 1997)). In a three-phase system, the requisite sensed signal contains two or more parameters that require a multi-sensing method. The tomographic technique for imaging a three-phase flow system is referred to as multi-modal tomography. Most of the tomographic techniques developed so far are for single-modal systems, which are not readily applicable to three-phase systems. There are three strategies to perform three-phase imaging using a tomographic technique: (1) by combination of two different single-modal sensing systems, (2) using an inherently multi-modal sensing system, and (3) by means of a single-modal sensing system having a reconstruction technique capable of differentiating between three phases in the object space (Warsito et al., 1999).
An example of the first approach is the use of electrical capacitance tomography combined with Gamma-ray tomography for imaging a multi-component composition of gas, oil and water in a pipeline (Johansen et al., 1996). Water has a permittivity approximately 40 times to that of oil, or 80 times to that of gas, so that it can be easily differentiated from gas and oil using the electrical capacitance tomography. The gas distribution is then determined using Gamma-tomography which is responsive to the large density difference between the gas and the liquids. George et al. (2000) used electrical resistance tomography (ERT) combined with Gamma densitometry tomography (GDT) to measure the gas and solid concentration (holdup) profiles in three-phase bubble column. They employed dry air as the gas phase, water with sodium nitrate as the liquid phase, and both polystyrene and glass-beads for the solid phase. Polystyrene is an electrical insulator like air but has a gamma attenuation coefficient similar to water, so that ERT is influenced by both the solid and gas phases while GDT is primarily sensitive to the gas alone. In the reconstruction, the gas volume fraction profile from GDT is subtracted from the insulating phase profile determined by ERT to get first-approximations of volume fractions of the three components. Glass has a gamma attenuation coefficient significantly different from both air and water, so that the first-order approximation does not hold in this case. Instead, the gamma attenuation formulas and the Maxwell-Hewitt conductivity relation for mixed materials must be solved as a set of simultaneous equations to reconstruct the volume fraction profiles of the three phases. However, the primary problem of using this approach arises because the measurements using the combined techniques are not conducted simultaneously in the same object domain, so that the reconstructed profiles may be severely distorted. In addition, this approach is complex to implement and impractical with respect to the overall cost.
The second approach uses a single sensing technique, which is an inherently multi-modal system capable of differentiating between two or more species in the object space. This approach is advantageous because all the information required is available using the same measuring technique and image reconstruction. An example of this approach is the dual-frequency ultrasonic method implemented by Warsito et al. (1995) to measure gas and solid concentration distributions in a three-phase bubble column. The gas and solid concentrations are obtained by measuring time-of-flights of ultrasonic pulse beams in high and low frequencies. In high frequency (10 MHzxcx9c), the response of the gas bubble on the time-of-flight is negligible, allowing the measurement of the solid phase alone. The gas phase is then calculated by subtracting the time-of-flight in lower frequency (1xcx9c3 MHz) by the time-of-flight due to the solid phase. Another technique using ultrasonic tomography with two-parameter sensing is implemented by Warsito et al. (1997, 1999) to measure the cross-sectional distributions of gas and solid concentrations in a gas-liquid-solid bubble column. They evaluated both the attenuation and the time-of-flight of a single ultrasonic pulse beam transmitted through the three-phase medium to construct a set of simultaneous equations involving gas and solid concentrations. The gas and the solid concentrations are obtained by solving the two equations. The major advantage of using this technique is that the gas and the solid concentrations are measured simultaneously. However, the technique is limited to relatively low concentrations of gas or solid, lower than 20% for both. The use of optical tomography coupled with spectral analysis and electrical tomography by exploring the potential of impedance spectroscopy and dielectric spectroscopy also shows potential for multi-modal tomography systems (Beck, 1995).
The third approach uses a single modal sensing system but with a reconstruction technique capable of differentiating between three species in the object domain. An example is the use of electrical capacitance tomography with a neural-network-based image reconstruction technique as proposed by Nooralahiyan and Hoyle (1997). They used a single layer feed forward neural network with double-step sigmoid function to replace a single-step sigmoid function in the neural-network computing to enable the identification of gas bubbles and water drops in oil environment. Though the network training is extremely time consuming, once well trained the neural network is very fast when used for prediction (reconstruction). However, the problem is that the feed forward neural network needs prior knowledge of the flow pattern for training before any measurements can be taken. This makes the technique impractical for real applications when training is not possible, particularly for the real time imaging of complex flows where the pattern is highly fluctuated and unknown before an exact image is obtained.
Principles of Electrical Capacitance Tomography (ECT)
An ECT system can be divided into three basic components: (1) capacitance sensor, (2) sensing electronics for data acquisition and (3) computer system for image reconstruction, interpretation and display. The capacitance sensor consists of a number of electrodes located on the periphery of the process vessel (FIG. 1). If the number of electrodes is ne, there will be ne(nexe2x88x921)/2 combinations of independent capacitance measurements between electrode-pairs. The measured capacitance is a function of the dielectric constant (permittivity) filling the space between the electrodes in the pair. Fundamental problems in ECT are associated with the process of image reconstruction from the measured capacitances based on Poisson""s equation given by (Xie et al., 1992, Xie, 1995)
∈(x,y)∇2"PHgr"(x,y)+∇∈(x,y)∇"PHgr"(x,y)=0xe2x80x83xe2x80x83(1)
where ∈(x,y) and "PHgr"(x,y) are the dielectric constant and the electrical potential distributions, respectively. This second-order, elliptic-type partial differential equation has either Dirichlet-type or Neumann-type typical boundary conditions. The dielectric field is given by
E=xe2x88x92∇"PHgr"(x,y)xe2x80x83xe2x80x83(2)
By applying the Gauss law, the charge at the detector electrode induced by the source electrode in the ith electrode-pair (see FIG. 1) can be calculated by the following expression:
Qi=∫∈(x,y)Exc2x7{circumflex over (n)}dlxe2x80x83xe2x80x83(3)
where xcex93i is a closed curve enclosing the detector electrode and {circumflex over (n)} is the unit normal vector to xcex93i. The capacitance between the electrode-pair, Ci, can be calculated by                               C          i                =                                            Q              i                                      Δ              ⁢                              xe2x80x83                            ⁢                              V                i                                              =                                    1                              Δ                ⁢                                  xe2x80x83                                ⁢                                  V                  i                                                      ⁢                                          ∮                                  Γ                  i                                            ⁢                              xe2x80x83                            ⁢                              ϵ                ⁢                                  xe2x80x83                                ⁢                                  (                                      x                    ,                    y                                    )                                ⁢                                  xe2x80x83                                ⁢                                  E                  ·                                      n                    ^                                                  ⁢                                  xe2x80x83                                ⁢                                  ⅆ                  l                                                                                        (        4        )            
where *Vi is the voltage difference between the source and the detector electrodes of electrode-pair i.
Equation (4) relates the dielectric constant (permittivity) distribution, ∈(x,y), to the measured capacitances Ci. That is, for a given medium distribution ∈(x,y) and the boundary conditions, the capacitances can be calculated (e.g., using a finite element method). This task is referred to as a forward problem. The image reconstruction process is an inverse problem involving the estimation of the permittivity distribution from the measured capacitances. However, as seen from this equation, it is difficult to obtain an explicit expression that relates the measured capacitances to the fraction and position of the dielectric components inside the pipe. Some reconstruction techniques developed for ECT are overviewed in the following section.
Image Reconstruction Techniques for ECT
Since there is no general solution method for the forward problem (in which the equation becomes non-linear), approximation methods are usually used. The most common one is using the so-called sensitivity model (Huang et al., 1989, Xie et al., 1992). The model is based on the electrical network superposition theorem in which the domain (the cross section of the sensor) is subdivided into a number of pixels, and the response of the sensor may be found as a sum (linear model) of interactions which will result when the permittivity of only one pixel in the domain is changed by a known amount. This is similar to the first order series expansion approach for xe2x80x98hard fieldxe2x80x99 tomography (Herman, 1980). However, this model only holds if the difference between the dielectric constants of the constituent materials being imaged is small (high-to-low permittivity ratio is less than 6). According to this model, eq. (4) can be written as a matrix expression as:
C=SGxe2x80x83xe2x80x83(5)
where C is the measured capacitance matrix, G is the image vector (permittivity distribution) and S is the so-called sensitivity map. Equation (5) is regarded as the linear forward projection (LFP). The sensitivity map is calculated as                               S          ij                =                              (                                          A                max                                            Δ                ⁢                                  xe2x80x83                                ⁢                                  ϵ                  0                                ⁢                                  A                  j                                                      )                    ⁢                      xe2x80x83                    ⁢                                                    C                i                                  (                                      H                    ,                    j                                    )                                            -                              C                i                L                                                                    C                i                H                            -                              C                i                L                                                                        (        6        )            
i=1, 2, . . . , M
j=1, 2, . . . , M
where Ci(H,j) is the capacitance value of the ith electrode-pair when pixel j inside the sensor is a material with high permittivity and the rest is material of low permittivity. M is the number of electrode-pair combinations, and N is the number of pixels in the object space. xcex94∈0 is the permittivity difference between the high and low permittivity materials. Axe2x80x2j is the area of element j and Amax is the largest element of all Aj. The sensitivity map is constructed by numerically solving equation (4) using a finite-element method (FEM) or by an experimental method using a dielectric rod placed at different positions inside the cross section of the sensor and measuring the capacitance between the electrode pair. The sensitivity map is necessary to solve the inverse problem using linear back projection (LBP) algorithm (Huang et al., 1989, Xie et al., 1992). Here, the image vector G is obtained from the following matrix transformation:
G=STCxe2x80x83xe2x80x83(7)
where ST is the matrix transpose of S.
The reconstructed image using LBP is blurred, demonstrating a smoothing effect on the sharp transitions between the different dielectric constants. This effect is typical in reconstruction results using the linear back-projection technique, because the gray level of each pixel in the image is formed from a set of overlapping projections. To reduce the blurring effect, a filter function as commonly used in filtered-back projection algorithm for linear tomography or a thresholding procedure is usually applied (Xie et al., 1992, Mwambela et al., 1997), however, only slight improvements have been achieved by these procedures.
A great demand for high quality images generated by ECT has led many researchers to develop reconstruction techniques based on iterative algorithms (Reinecke and Mewes, 1996, 1997, Isaken, 1996, Yang, 1997, Liu et al., 1997, Su et al., 2000). Iterative algorithms may be classified into two groups: an algebraic reconstruction technique (ART) type and an optimization type. The ART type is a modification of the original technique developed for linear tomography by combining the linear reconstruction technique with the sensitivity model. The image reconstruction problem with using ART is estimating the image vector (permittivity distribution) G such that the estimated capacitances C*(G)xe2x89xa6C, given a measurement vector C. The following form generally represents ART type algorithms:
G(k+1)=G(k)+xcex1(k)ST(Cxe2x88x92C*(G(k)))xe2x80x83xe2x80x83(8)
where G(k) is the estimated permittivity vector in the k-th iteration, xcex1 is a relaxation factor, also called a gain factor or weighting factor. The initial estimate of the permittivity vector G(0) is calculated using LBP by eq. (7). For the sake of convenience, eq. (8) is normalized as
{tilde over (G)}(k+1)={tilde over (G)}(k)+xcex1(k){tilde over (S)}T({tilde over (C)}xe2x88x92{tilde over (C)}*({tilde over (G)}(k)))xe2x80x83xe2x80x83(9)
where                                                         G              ~                        j                    =                                                    G                j                            -                              ϵ                                  0                  ⁢                  L                                                                                    ϵ                                  0                  ⁢                  H                                            -                              ϵ                                  0                  ⁢                  L                                                                    ,                  xe2x80x83                ⁢                                            S              ~                                      i              ⁢                              xe2x80x83                            ⁢              j                                =                                    S                              i                ⁢                                  xe2x80x83                                ⁢                j                                                                    ∑                j                            ⁢                              S                                  i                  ⁢                                      xe2x80x83                                    ⁢                  j                                                                    ,                  
                ⁢                                            S              ~                                      i              ⁢                              xe2x80x83                            ⁢              j                        T                    =                                                    S                ~                                            j                ⁢                                  xe2x80x83                                ⁢                i                                      =                                          S                                  j                  ⁢                                      xe2x80x83                                    ⁢                  i                                                                              ∑                  i                                ⁢                                  S                                      i                    ⁢                                          xe2x80x83                                        ⁢                    j                                                                                      ,                  xe2x80x83                ⁢                                            C              ~                        i                    =                                                    C                i                            -                              C                i                L                                                                    C                i                H                            -                              C                i                L                                                                        (        10        )            
i=1,2, . . . , M
j=1,2, . . . , N
For simplicity, in the rest of this paper the normalization notation ({tilde over ( )}) is omitted, and all the above parameters, otherwise stated, use the normalized forms.
Many ART type algorithms for ECT used by researchers differ slightly from each other in the structures of the relaxation factor xcex1 and the procedure used to correct the estimated permittivity values. Commercially available reconstruction software developed by PTL (UK, 1999) uses a constant relaxation factor and LFP technique (eq. 5) to solve the forward problem in the correction term for the estimated permittivity, i.e.:
C*(G(k))=SG(k).xe2x80x83xe2x80x83(11)
This is simply LBP used in an iterative manner, similar to iterative back projection techniques for linear tomography. Reinecke and Mewes (1996, 1997) use a correction factor       (                  ∑        j            ⁢              S                  i          ⁢                      xe2x80x83                    ⁢          j                2              )        -    1  
instead of the relaxation factor, similar to ART for linear tomography (Gordon et al., 1970), and use a finite difference technique to solve the forward problem in the correction term to replace the LFP method. It should be noted that, unlike Gordon ART for xe2x80x98hard fieldxe2x80x99 (linear) tomography, in electrical capacitance tomography problems the electrical field spreads all over the measured cross-section allowing no distinct border between the pixels which contribute to the ith pair capacitance measurement and those which do not, due to the xe2x80x98soft fieldxe2x80x99 effect. Thus, for ECT, all of the capacitance data is typically used to update the pixel values. In this regard, the technique can also be referred to as simultaneous iterative reconstruction technique (SIRT) as used by Su et al. (2000). The technique used by Reinecke and Mewes (1997) has improved the accuracy of the reconstructed result, though the finite difference technique significantly increases the computation time. To improve the convergence performance, Su et al. (2000) determined a priori the relaxation factor in the form of:                                           α                          (              k              )                                =                                    α              0                        +                          β              k                                      ,                            (        12        )            
where xcex10 and xcex2 are positive constants. They also used LFP to solve the forward problem in the correction term in eq. (10). The choice of the relaxation factor has significantly improved the convergence performance, though there is still a smoothing effect in the region between different permittivities. Liu et al. (1999) used an optimal relaxation factor calculated by minimizing the error function e(k+1)=∥Cxe2x88x92SG(k+1)∥2 so that             ∂              e                  (                      k            +            1                    )                            ∂              α                  (          k          )                      =  0.
Yang et al. (1999) determined the gain factor based on suitable convergence criterion:
xcex1(k)=2/xcexmax,xe2x80x83xe2x80x83(13)
where xcexmax is the maximum eigenvalue of STS. They replaced LFP with FEM to solve the forward problem in the correction term to improve both the image quality and the convergence rate, though this significantly increases the computation load.
The second group of iterative techniques is based on an optimization algorithm in which one or more objective functions are minimized or maximized. The objective functions, typically the least squared error sum and/or the maximum entropy function, are measures of the xe2x80x98likelihoodxe2x80x99 with respect to the definition of the problem. Isaken and Nordtvetd (1996) proposed an optimization technique for ECT called model-based reconstruction (MOR). The MOR technique minimizes the least squared error sum between the measured capacitance and the estimation value calculated using model-based parameter chosen to represent the image pattern. Thus, a prior knowledge of the image pattern (permittivity distribution) is required for the pattern representation before the reconstruction is done. Mwambela et al. (1997) used the maximum entropy function to generate threshold values to address the xe2x80x98smoothingxe2x80x99 effect on LBP results.
With regard to minimizing an objective function like MOR, ART is a kind of optimization technique aimed at finding a least squared error between the measured and the estimated capacitances. Such criteria, however, does not necessarily determine the image vector G because image reconstruction is an ill-posed problem as there are fewer independent measurements than unknown pixel values. Therefore, there may be more than one possible alternative image. Another reason why such a solution is not necessarily very good is that the least squared criterion does not contain any information regarding the nature of a xe2x80x98desirablexe2x80x99 solution (Herman, 1980). Therefore, more than one objective function is required to be minimized simultaneously to choose the xe2x80x98best compromise-solutionxe2x80x99 from possible alternatives.
In this work, a multi-criteria optimization based image reconstruction technique for imaging gas-liquid two-phase flows using electrical capacitance tomography is developed. The reconstruction technique is a combination of multi-criteria optimization image reconstruction technique (Wang et al., 1996, Wang, 1997, 1998) for linear tomography and the LBP technique. The multi-criteria optimization image reconstruction problem is then solved using modified Hopfield model dynamic neural-network computing. It is believed that this is the first application of multi-criteria optimization technique to a xe2x80x98soft fieldxe2x80x99 tomography problem.
The reconstruction technique for the two-phase system is extended to three-phase system by introducing a double-step sigmoid in the modified Hopfield network. By using the double-step sigmoid function, the calculation is permitted to converge to three-different stable regions in the output space corresponding to the three-different phases, enabling differentiation among the three components. This technique has a major advantage over the feed-forward neural network based image reconstruction technique as used by Nooralahiyan and Hoyle (1997) because there is no network training needed. Thus, it can be applied to any multiphase system without prior knowledge on the flow pattern.
The image reconstruction technique of the present invention is based upon multi-criteria optimization using an analog neural network (NN-MOIRT). This reconstruction technique is a combination of multi-criteria optimization image reconstruction techniques for linear tomography, and the so-called linear back projection (LBP) technique commonly used for capacitance tomography. The multi-criteria optimization image reconstruction problem is solved using Hopfield model dynamics neural-network computing.
A plurality of sensors is placed around the periphery of a conduit through which flows a two- or three-phase fluid system. The sensors collect image measurements of the fluid between sensors. That is to say, the plurality of sensors will provide n(nxe2x88x921)/2 image measurements, where n is the number of sensors. The image measurements shall collectively be referred to as xe2x80x98dataxe2x80x99 or xe2x80x98image measurement dataxe2x80x99.
The present invention presents a novel method for image reconstruction comprising obtaining data, processing the data by application of a linear back projection algorithm with a Hopfield neural network, and displaying the processed data on a display device.
Suitable forms of data include, but are not limited to, capacitance data, conductance data, x-ray data, gamma-ray data, ultrasonic data and optical data. It is preferred that the data is capacitance data or conductance data. It is most preferred that the data is capacitance data.
The linear back projection algorithm may be an iterative back projection algorithm (ILBP). The linear back projection algorithm preferably comprises at least one function to minimize. It is most preferred that the linear back projection algorithm comprises three functions to minimize: a negative entropy function, a weighted square error function and a sum of non-uniformity and peakedness functions.
It is preferred that the Hopfield model neural network is modified to include a penalty function that permits a temporary increase in the descending evolution of the Hopfield network energy to allow escape from local minima.
Suitable display devices include, but are not limited to: monitors, projectors, printers and plotters. Other display devices may be suitable for displaying the data and will be known to those skilled in the art.
The present invention provides a method for obtaining a cross-sectional image of a two-phase fluid flowing through a conduit. A sensor is comprised of a transmitter and a receiver. The method comprises acquiring image measurement data from at least two sensors peripherally located on the interior of the conduit, the at least two sensors defining a cross-section of the conduit. The image measurement data is next sent to a processing unit. The image measurement data processing comprises the simultaneous minimization of a negative entropy function of an image given by             f      1        ⁢          (      G      )        =            γ      1        ⁢                  ∑                  j          =          1                N            ⁢                        G          j                ⁢                  xe2x80x83                ⁢        ln        ⁢                  xe2x80x83                ⁢                  G          j                    
where xcex31 is a normalized constant between 0 and 1; a weighted square error function given by             f      2        ⁢          (      G      )        =            1      2        ⁢          xe2x80x83        ⁢          γ      2        ⁢                  ∑                  i          =          1                M            ⁢                        (                                                    ∑                                  j                  =                  1                                N                            ⁢                                                S                                      i                    ⁢                                          xe2x80x83                                        ⁢                    j                                                  ⁢                                  G                  j                                                      ⁢                          xe2x80x83                        -                          C              i                                )                2            
where xcex32 is a normalized constant between 0 and 1; and a sum of non-uniformity and peakedness functions given by             f      3        ⁢          (      G      )        =            1      2        ⁢          xe2x80x83        ⁢                  γ        3            ⁢              (                                            G              T                        ⁢            X            ⁢                          xe2x80x83                        ⁢            G                    +                                    G              T                        ⁢            G                          )            
where xcex33 is a normalized constant between 0 and 1. The simultaneous minimization is accomplished by application of a Hopfield model neural network given by       E    ⁢          (      G      )        =                    ∑                  l          =          1                3            ⁢                        ω          l                ⁢                              f            l                    ⁢                      (            G            )                                +                  ∑                  i          =          1                M            ⁢              Ψ        ⁢                  xe2x80x83                ⁢                  (                      z            i                    )                      +                  ∑                  j          =          1                N            ⁢                        1                      R            j                          ⁢                              ∫            0                          G              j                                ⁢                                                    f                Σ                                  -                  1                                            ⁢                              (                G                )                                      ⁢                          ⅆ              G                                          
wherein normalized permittivity values are mapped into output variables using a single sigmoid function given by                     f        Σ            ⁢              (                  u          j                )              =          1              1        +                  exp          ⁢                      xe2x80x83                    ⁢                      (                                          -                β                            ⁢                              xe2x80x83                            ⁢                              u                j                                      )                                ,
to force convergence to a binary output. The Hopfield model neural network further comprising a penalty function xcexa8(zi) defined as       Ψ    ⁢          xe2x80x83        ⁢          (              z        i            )        =      {                            0                                                    if              ⁢                              xe2x80x83                            ⁢                              z                i                                      ≤            0                                                            α            ⁢                          xe2x80x83                        ⁢                          z              i                                                                          if              ⁢                              xe2x80x83                            ⁢                              z                i                                       greater than             0                              
where       z    i    =                    ∑                  j          =          1                N            ⁢                        S                      i            ⁢                          xe2x80x83                        ⁢            j                          ⁢                  G          j                      ⁢          xe2x80x83        -          C      i      
with xcex1(t)=xcex10+"xgr" exp(xe2x88x92xcex7t) where xcex10, "xgr", and xcex7 are positive constants. The penalty function prevents convergence at local minima. The output variables are then sent to a display device to display a cross-sectional image of the two-phase fluid flowing through the conduit.
The two-phase fluid flowing through the interior of the conduit may consist of any combination of solids and gases; solids and liquids; liquids and gases; two different solids; two different liquids; or two different gases.
The present invention also provides for a method for obtaining a cross-sectional image of a three-phase fluid flowing through a conduit having an interior. The method comprises acquiring image measurement data from at least two sensors peripherally located on the interior of the conduit, the at least two sensors defining a cross-section of said conduit. The image measurement data is then sent to a processing unit. The image measurement data processing comprises the simultaneous minimization of a negative entropy function of an image given by             f      1        ⁡          (      G      )        =            γ      1        ⁢                  ∑                  j          =          1                N            ⁢                        G          j                ⁢                  xe2x80x83                ⁢        ln        ⁢                  xe2x80x83                ⁢                  G          j                    
where xcex31 is a normalized constant between 0 and 1; a weighted square error function given by             f      2        ⁡          (      G      )        =            1      2        ⁢          xe2x80x83        ⁢          γ      2        ⁢                  ∑                  i          =          1                M            ⁢                        (                                                    ∑                                  j                  =                  1                                N                            ⁢                                                S                                      i                    ⁢                                          xe2x80x83                                        ⁢                    j                                                  ⁢                                  G                  j                                                      ⁢                          xe2x80x83                        -                          C              i                                )                2            
where xcex32 is a normalized constant between 0 and 1; and a sum of non-uniformity and peakedness functions given by             f      3        ⁡          (      G      )        =            1      2        ⁢          xe2x80x83        ⁢                  γ        3            ⁡              (                                            G              T                        ⁢            X            ⁢                          xe2x80x83                        ⁢            G                    +                                    G              T                        ⁢            G                          )            
where xcex33 is a normalized constant between 0 and 1. The simultaneous minimization is accomplished by application of a Hopfield model neural network given by       E    ⁡          (      G      )        =                    ∑                  l          =          1                3            ⁢                        ω          l                ⁢                              f            l                    ⁡                      (            G            )                                +                  ∑                  i          =          1                M            ⁢              Ψ        ⁢                  xe2x80x83                ⁢                  (                      z            i                    )                      +                  ∑                  j          =          1                N            ⁢                        1                      R            j                          ⁢                              ∫            0                          G              j                                ⁢                                                    f                Σ                                  -                  1                                            ⁡                              (                G                )                                      ⁢                          ⅆ              G                                          
wherein normalized permittivity values are mapped into output variables using a double sigmoid function given by                     f        Σ            ⁡              (                  u          j                )              =                            f                      Σ            ⁢            1                          ⁡                  (                      u            j                    )                    +                        f                      Σ            ⁢            2                          ⁡                  (                      u            j                    )                      ,
where                     f                  Σ          ⁢          1                    ⁡              (                  u          j                )              =                  ρ        1                    1        +                  exp          ⁢                      xe2x80x83                    ⁢                      (                                          -                                                      β                    ⁢                                          xe2x80x83                                                        1                                            ⁢                              (                                                      u                    j                                    +                                      ξ                    1                                                  )                                      )                                ,                    f                  Σ          ⁢          2                    ⁡              (                  u          j                )              =                  ρ        2                    1        +                  exp          ⁢                      xe2x80x83                    ⁢                      (                                          -                                                      β                    ⁢                                          xe2x80x83                                                        2                                            ⁢                              (                                                      u                    j                                    +                                      ξ                    2                                                  )                                      )                                ,
xcfx811 and xcfx812 are positive constants, xcex21 and xcex22 are steepness gains of sigmoid functions ƒxcexa31 and ƒxcexa32, and "xgr"1 and "xgr"2 are parameters, where the double sigmoid function purpose is to force convergence to a tertiary output. The Hopfield model further comprises a penalty function xcexa8(zi) defined as       Ψ    ⁢          xe2x80x83        ⁢          (              z        i            )        =      {                            0                                                    if              ⁢                              xe2x80x83                            ⁢                              z                i                                      ≤            0                                                            α            ⁢                          xe2x80x83                        ⁢                          z              i                                                                          if              ⁢                              xe2x80x83                            ⁢                              z                i                                       greater than             0                              
where       z    i    =                    ∑                  j          =          1                N            ⁢                        S                      i            ⁢                          xe2x80x83                        ⁢            j                          ⁢                  G          j                      ⁢          xe2x80x83        -          C      i      
with xcex1(t)=xcex10+"xgr" exp(xe2x88x92xcex7t) where xcex10, "xgr", and xcex7 are positive constants. The penalty function prevents convergence at local minima. The output variables are then sent to a display device to display a cross-sectional image of the three-phase fluid flowing through the conduit.
The three-phase fluid flow consists of any combination of solids, liquids and/or gases.
The present invention further provides a method for obtaining a cross-sectional image of a multiphase fluid flowing through a conduit having an interior using a Hopfield type neural network by acquiring image measurement data for fluid flowing through the conduit. The image measurement data is acquired by at least two sensors peripherally located on the interior of the conduit, the at least two sensors defining a cross-section of the conduit. The Hopfield model neural network is initialized by (1) choosing an initial state of neurons; (2) setting at least one steepness gain factor; (3) initializing a penalty parameter; and (4) initializing an initial gain factor. The image vector is updated accordingly. The iterative process stops when |Gj(t+xcex94t)xe2x88x92Gj(t)|2xe2x89xa6e for all neurons, where |Gj(t+xcex94t)xe2x88x92Gj(t)|2 is error and e is a termination scalar set by an operator. The image vector is then displayed on a suitable display device.
It is preferred that the initial state of neurons is set to zero. For a two-phase fluid flow, it is preferred that the steepness gain factor is set to two. For a three-phase fluid flow, it is preferred that both steepness gain factors are set to 2. It is further preferred that the initial state of the gain factor and the initial state of the penalty factor are each experimentally determined.
For a two-phase fluid flow, it is preferred that the image vector is updated by calculating coefficients of objective functions:             γ      1              (                  t          +                      Δ            ⁢                          xe2x80x83                        ⁢            t                          )              =                  [                              ∑                          j              =              1                        N                    ⁢                                                    G                j                            ⁡                              (                t                )                                      ⁢                          xe2x80x83                        ⁢            ln            ⁢                          xe2x80x83                        ⁢                                          G                j                            ⁡                              (                t                )                                                    ]                    -        1              ,            γ      2              (                  t          +                      Δ            ⁢                          xe2x80x83                        ⁢            t                          )              =                  [                              1            2                    ⁢                                    "LeftDoubleBracketingBar"                                                SG                  ⁡                                      (                    t                    )                                                  -                C                            "RightDoubleBracketingBar"                        2                          ]                    -        1              ,      
    ⁢                    γ        3                  (                      t            +                          Δ              ⁢                              xe2x80x83                            ⁢              t                                )                    =                        [                                                    1                2                            ⁢                                                G                  T                                ⁡                                  (                  t                  )                                            ⁢              s              ⁢                              xe2x80x83                            ⁢                              G                ⁡                                  (                  t                  )                                                      +                                          1                2                            ⁢                                                G                  T                                ⁡                                  (                  t                  )                                            ⁢                              G                ⁡                                  (                  t                  )                                                              ]                          -          1                      ;  
updating weights xcfx891,xcfx892,xcfx893 for every iteration as follows:             ω      i              (                  t          +                      Δ            ⁢                          xe2x80x83                        ⁢            t                          )              =                            Δω          1                      (            t            )                          /                  Δω          i                      (            t            )                                                            ∑            3                                i            =            1                          ⁢                              Δω            1                          (              t              )                                /                      Δω            i                          (              t              )                                            ,
xcex94xcfx891(t)=ƒi(G(t+xcex94t))xe2x88x92ƒi(G(t)), (i=1,2,3) where f1 is a negative entropy function, f2 is a weighted square error function and f3 is a sum of non-uniformity and peakedness functions of an image; and calculating an image vector Gj(t+xcex94t), where Gj(t+xcex94t)=Gj(t)+ƒxe2x80x2xcexa3(u)uxe2x80x2j(t)xcex94t, ƒxe2x80x2xcexa3(u)=dƒxcexa3(uj)/duj and uxe2x80x2j(t)=duj(t)/dt.
For a three-phase fluid flow, it is preferred that the image vector is updated by calculating coefficients of objective functions:             γ      1              (                  t          +                      Δ            ⁢                          xe2x80x83                        ⁢            t                          )              =                  [                                            ∑              N                                      j              =              1                                ⁢                                                    G                j                            ⁢                              (                t                )                                      ⁢            ln            ⁢                          xe2x80x83                        ⁢                                          G                j                            ⁢                              (                t                )                                                    ]                    -        1              ,            γ      2              (                  t          +                      Δ            ⁢                          xe2x80x83                        ⁢            t                          )              =                  [                              1            2                    ⁢                                    "LeftDoubleBracketingBar"                                                SG                  ⁢                                      (                    t                    )                                                  -                C                            "RightDoubleBracketingBar"                        2                          ]                    -        1              ,      
    ⁢                    γ        3                  (                      t            +                          Δ              ⁢                              xe2x80x83                            ⁢              t                                )                    =                        [                                                    1                2                            ⁢                                                G                  T                                ⁢                                  (                  t                  )                                            ⁢                              sG                ⁢                                  (                  t                  )                                                      +                                          1                2                            ⁢                                                G                  T                                ⁢                                  (                  t                  )                                            ⁢                              G                ⁢                                  (                  t                  )                                                              ]                          -          1                      ;  
updating weights xcfx891,xcfx892,xcfx893 for every iteration as follows:             ω      i              (                  t          +                      Δ            ⁢                          xe2x80x83                        ⁢            t                          )              =                            Δω          1                      (            t            )                          /                  Δω          i                      (            t            )                                                            ∑            3                                i            =            1                          ⁢                              Δω            1                          (              t              )                                /                      Δω            i                          (              t              )                                            ,
xcex94xcfx89i(t)=ƒi(G(t+xcex94t))xe2x88x92ƒi(G(t)), (i=1,2,3) where f1 is a negative entropy function, f2 is a weighted square error function and f3 is a sum of non-uniformity and peakedness functions of an image; and calculating an image vector Gj(t+xcex94t), where Gj(t+xcex94t)=ƒxcexa3(uj(t+xcex94t))=Gj(t)+[ƒxcexa31(uj)(1xe2x88x92ƒxcexa31(uj))+ƒxcexa32(uj)(1xe2x88x92ƒxcexa32(uj))]uxe2x80x2j(t)xcex94t,                     f                  Σ          ⁢                      xe2x80x83                    ⁢          1                    ⁡              (                  u          j                )              =                  ρ        1                    1        +                  exp          ⁡                      (                          -                                                β                  1                                ⁡                                  (                                                            u                      j                                        +                                          ξ                      1                                                        )                                                      )                                ,      
    ⁢                    f                  Σ          ⁢                      xe2x80x83                    ⁢          2                    ⁡              (                  u          j                )              =                  ρ        2                    1        +                  exp          ⁡                      (                          -                                                β                  2                                ⁡                                  (                                                            u                      j                                        +                                          ξ                      2                                                        )                                                      )                                ,
xcfx811={tilde over (∈)}0M, and xcfx812=1xe2x88x92{tilde over (∈)}0M, where {tilde over (∈)}0M is the normalized relative permittivity of the medium-permittivity material.